Solutions: Mathematica 6 Lab 5

Problem 1: Invaded Planetary System

Show[Graphics[{
{AbsolutePointSize[10],Point[{0,0}]},
{Red,Circle[{0,0},1]},
{Red,AbsolutePointSize[4],Point[{1,0}]},
{Blue,Line[{{-10,4},{10,4}}]},
{Blue,AbsolutePointSize[8],Point[{10,4}]}
}],AspectRatio->Automatic]

UnperturbedOrbit=NDSolve[{x''[t]==-((9x[t])/(((((x[t])^2)+((y[t])^2)))^(3/2))),y''[t]==-((9y[t])/(((((x[t])^2)+((y[t])^2)))^(3/2))),x[0]==1,y[0]==0,x'[0]==0,y'[0]==3},{x[t],y[t]},{t,0,3}]

ShowUnperturbed=ParametricPlot[Evaluate[{x[t],y[t]}/.UnperturbedOrbit],{t,0,2},
PlotRange->All,
 AspectRatio->Automatic,
PlotStyle->Red]

PerturbedOrbit=NDSolve[{x''[t]==-((9x[t])/(((((x[t])^2)+((y[t])^2)))^(3/2)))+((1/2(10-x[t]))/((((((10-x[t]))^2)+(((4-y[t]))^2)))^(3/2))),
y''[t]==-((9y[t])/(((((x[t])^2)+((y[t])^2)))^(3/2)))+((1/2(4-y[t]))/((((((10-x[t]))^2)+(((4-y[t]))^2)))^(3/2))),x[0]==1,y[0]==0,x'[0]==0,y'[0]==3},{x[t],y[t]},{t,0,80}]

ShowPerturbed=ParametricPlot[Evaluate[{x[t],y[t]}/.PerturbedOrbit],{t,0,80},
PlotRange->All,
 AspectRatio->Automatic,
PlotStyle->Red]

ImpactedOrbit=NDSolve[{x''[t]==-((9x[t])/(((((x[t])^2)+((y[t])^2)))^(3/2)))+((1/2(10-(t/4)-x[t]))/((((((10-(t/4)-x[t]))^2)+(((4-y[t]))^2)))^(3/2))),
y''[t]==-((9y[t])/(((((x[t])^2)+((y[t])^2)))^(3/2)))+((1/2(4-y[t]))/((((((10-(t/4)-x[t]))^2)+(((4-y[t]))^2)))^(3/2))),x[0]==1,y[0]==0,x'[0]==0,y'[0]==3},{x[t],y[t]},{t,0,80}]

ShowImpacted=ParametricPlot[Evaluate[{x[t],y[t]}/.ImpactedOrbit],{t,0,80},
PlotRange->All,
 AspectRatio->Automatic,
PlotStyle->Red]

ReverseImpact=NDSolve[{x''[t]==-((9x[t])/(((((x[t])^2)+((y[t])^2)))^(3/2)))+((1/2(-10+(t/4)-x[t]))/((((((-10+(t/4)-x[t]))^2)+(((4-y[t]))^2)))^(3/2))),
y''[t]==-((9y[t])/(((((x[t])^2)+((y[t])^2)))^(3/2)))+((1/2(4-y[t]))/((((((-10+(t/4)-x[t]))^2)+(((4-y[t]))^2)))^(3/2))),x[0]==1,y[0]==0,x'[0]==0,y'[0]==3},{x[t],y[t]},{t,0,80}]

ShowReverse=ParametricPlot[Evaluate[{x[t],y[t]}/.ReverseImpact],{t,0,80},
PlotRange->All,
 AspectRatio->Automatic,
PlotStyle->Red]

Problem 2:

Hyperlink["Van der Pol","http://www.scholarpedia.org/article/Van_der_Pol_Oscillator"]

Clear[b,x,μ,t]

b[x_,μ_]:=μ((x^2)-1)

solv1=NDSolve[{x''[t]+b[x[t],1/8]x'[t]+x[t]==0,x[0]==1,x'[0]==0},x[t],{t,0,30}]

graph1=Plot[x[t]/.solv1,{t,0,30}]

solv2=NDSolve[{x''[t]+b[x[t],1/8]x'[t]+x[t]==0,x[0]==4,x'[0]==0},x[t],{t,0,30}]

graph2=Plot[x[t]/.solv2,{t,0,30}]

solv3=NDSolve[{x''[t]+b[x[t],8]x'[t]+x[t]==0,x[0]==1,x'[0]==0},x[t],{t,0,30}]

graph3=Plot[x[t]/.solv3,{t,0,30}]

phase1=NDSolve[{x'[t]==y[t],y'[t]+8(((x[t])^2)-1)y[t]+x[t]==0,x[0]==0,y[0]==1},{x[t],y[t]},{t,0,30}]

phasemap1=ParametricPlot[{x[t],y[t]}/.phase1,{t,0,30},
PlotRange->All,
 AspectRatio->Automatic,
PlotStyle->Red,Ticks->None]

phase2=NDSolve[{x'[t]==y[t],y'[t]+8(((x[t])^2)-1)y[t]+x[t]==0,x[0]==0,y[0]==12},{x[t],y[t]},{t,0,30}]

phasemap2=ParametricPlot[{x[t],y[t]}/.phase2,{t,0,30},
PlotRange->All,
 AspectRatio->Automatic,PlotStyle->Black,Ticks->None]

phase3=NDSolve[{x'[t]==y[t],y'[t]+8(((x[t])^2)-1)y[t]+x[t]==0,x[0]==2.1,y[0]==17.2},{x[t],y[t]},{t,0,30}]

phasemap3=ParametricPlot[{x[t],y[t]}/.phase3,{t,0,30},
PlotRange->All,
 AspectRatio->Automatic,
PlotStyle->Blue,Ticks->None]

Show[phasemap1, phasemap2, phasemap3]

Problem 3: Oscillators, Coupled in Various Ways

Uncoupled=NDSolve[{x''[t]+(2^2)x[t]==0,y''[t]+(5^2)y[t]==0,x[0]==1,y[0]==0,x'[0]==0,y'[0]==2},{x[t],y[t]},{t,0,100}]

ParametricPlot[{x[t],y[t]}/.Uncoupled,{t,0,100}, AspectRatio->Automatic,PlotRange->{{-1,1},{-.4,.4}}]

DampedUncoupled=NDSolve[{x''[t]+(2^2)x[t]==-.1x'[t],y''[t]+(5^2)y[t]==-.1y'[t],x[0]==1,y[0]==0,x'[0]==0,y'[0]==2},{x[t],y[t]},{t,0,100}]

ParametricPlot[{x[t],y[t]}/.DampedUncoupled,{t,0,100}, AspectRatio->Automatic,PlotRange->{{-1,1},{-.4,.4}}]

UndampedCoupled=NDSolve[{x''[t]+(2^2)x[t]==5y[t],y''[t]+(5^2)y[t]==5x[t],x[0]==1,y[0]==0,x'[0]==0,y'[0]==2},{x[t],y[t]},{t,0,100}]

ParametricPlot[{x[t],y[t]}/.UndampedCoupled,{t,0,100}, AspectRatio->Automatic]

UndampedCoupled2=NDSolve[{x''[t]+(2^2)x[t]==9y[t],y''[t]+(5^2)y[t]==9x[t],x[0]==1,y[0]==0,x'[0]==0,y'[0]==2},{x[t],y[t]},{t,0,100}]

ParametricPlot[{x[t],y[t]}/.UndampedCoupled2,{t,0,100}, AspectRatio->Automatic]

UndampedCoupled3=NDSolve[{x''[t]+(2^2)x[t]==9.999y[t],y''[t]+(5^2)y[t]==9.999x[t],x[0]==1,y[0]==0,x'[0]==0,y'[0]==2},{x[t],y[t]},{t,0,100}]

ParametricPlot[{x[t],y[t]}/.UndampedCoupled3,{t,0,100}, AspectRatio->Automatic]

DampedCoupled=NDSolve[{x''[t]+(2^2)x[t]==-.005x'[t]+9y[t],y''[t]+(5^2)y[t]==-.005y'[t]9+x[t],x[0]==1,y[0]==0,x'[0]==0,y'[0]==2},{x[t],y[t]},{t,0,100}]

ParametricPlot[{x[t],y[t]}/.DampedCoupled,{t,0,100}, AspectRatio->Automatic, PlotRange->{{-1.3,1.3},{-.5,.5}}]

GyroCoupled=NDSolve[{x''[t]+2y'[t]+(2^2)x[t]==0,y''[t]-2x'[t]+(5^2)y[t]==0,x[0]==1,y[0]==0,x'[0]==0,y'[0]==2},{x[t],y[t]},{t,0,100}]

ParametricPlot[{x[t],y[t]}/.GyroCoupled,{t,0,100}, AspectRatio->Automatic]

HybridCoupled=NDSolve[{x''[t]+.6y'[t]+(2^2)x[t]==0,y''[t]-.6x[t]+(5^2)y[t]==0,x[0]==1,y[0]==0,x'[0]==0,y'[0]==2},{x[t],y[t]},{t,0,100}]

ParametricPlot[{x[t],y[t]}/.HybridCoupled,{t,0,100}, AspectRatio->Automatic]